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A key observable signature of integrability—of the existence of infinitely many “higher” conservation 
laws Qj] —in a system supporting solitons is the fact that a collision between solitons does not change their 
shape or size. But then, if solitons meet on top of a strong integrahility-breaking barrier, one would expect the 
solitons to undergo some process consistent with energy conservation but not with higher conservation laws, 
such as the larger soliton cannibalizing the smaller one. However, here we show that when a strongly-coupled 
“breather” 0] of the integrahle nonlinear Schrddinger equation is scattered off a strong barrier, the solitons 
constituting the breather separate but survive the collision: as we launch a breather with a fixed impact speed 
at barriers of lower and lower height, at first all constituent solitons are fully reflected; then, at a critical barrier 
height, the smallest soliton gets to be fully transmitted, while the other ones are still fully reflected. This persists 
as the barrier is lowered some more until, at another critical height, the second smallest soliton begins to be 
fully transmitted as well, etc., resulting in a staircase-like transmission plot, with quantized plateaus. We show 
how this effect makes tangible the inverse scattering transform: the powerful, but otherwise physically opaque 
mathematical formalism for solving completely integrable partial differential equations. Furthermore, if such 
collisions can be realized experimentally in ultracold Bose gases sa , they could become a basis of improved 
atom interferometers mi]. 


Completely integrable partial differential equations (PDEs) 
have played a central conceptual role in studies of nonlin¬ 
ear dynamical systems, such as those of integrability-to-chaos 
transition and thermalization H. At the same time, many 
physical systems are well-modeled by them, so that if a PDE 
supports solitons, they are also objects of interest in their 
own right. Here we focus on the ID nonlinear Schrodinger 
equation (NLSE), Eq. (IS2I) below. As one of universal equa¬ 
tions describing envelope dynamics of quasi-monochromatic 
waves in weakly nonlinear dispersive media with negligible 
dissipation iQl, the NLSE describes a variety of systems, in¬ 
cluding light propagating in nonlinear planar waveguides and 
optical fibers H; small-amplitude gravity waves on the sur¬ 
face of deep inviscid water ||3l; the Langmuir waves in hot 
plasmas iQl ; storage and transfer of vibrational energy in a- 
helix proteins via Davydov’s solitons 1^; and, of special rel¬ 
evance to us, ultracold quantum Bose gases in the mean-field 
regime in highly elongated traps that make them effectively 
one-dimensional (ID). In particular, since their first ex¬ 
perimental realizations bright matter-wave solitons 

have also been studied in the context of matter-wave inter¬ 
ferometry. Their use may improve sensitivity by several or¬ 
ders of magnitude fiH], especially in precise force sensing 
jg H^ and in the measurement of small magnetic field gradi¬ 
ents ||a[I2]- Much-studied are the soliton collisions with po¬ 
tential barriers, resulting i n sp litting ifTsl - l^ and subsequent 
recombination of solitons . 

Here we will be interested in the regimes where a single 
soliton does not measurably split or radiate its norm away due 
to the collision; when it is, for all intents and purposes, ei¬ 
ther completely reflected or completely transmitted, in its en¬ 
tirety. Since we consider the classical held description, which 
is mean-field theory in the matter-wave context, this corre¬ 
sponds to the following regime: let the velocity of the soliton 
center of mass be u; let be the total number of particles. 
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each of mass m, in the soliton; and let i?rest/-^a be the soliton 
energy per particle in the frame in which its center of mass is 
at rest (so that this is also the interatomic interaction energy 
per particle). The regime of interest is defined by the center- 
of-mass kinetic energy per particle, i being 

much less than L^resiZ-^a- The solitary wave is then energet¬ 
ically allowed to lose a certain number of particles, AN, but 
the fraction it is allowed to lose, AN/N^, goes down as 1/N^ 
if we increase the N^. while keeping v constant (see the Sup- 
plemetary materials). The plot of transmission vs. incident 
kinetic energy will be almost a step function, with quantized 
plateaus at 0 and 1. 

But now suppose we scatter a composite soliton, consist¬ 
ing of, e.g., two solitons propagating together, one on top of 
the other. The mental picture outlined above would suggest 
that during the collision, the larger constituent soliton (i.e. its 
solitary-wave descendant in the nonintegrable regime) should 
absorb at least some portion of the particles from the smaller 
one. Not only is this process energetically favorable, but since 
the interactions within the constituent solitons are of the same 
order as the interactions between the constituent solitons, one 
would expect such processes to occur as soon as there are no 
conservation laws to forbid them. 

Nevertheless, we will now show that there are collision 
regimes when the constituent solitons are in fact preserved', 
and this despite the fact that the collision process is strongly 
nonperturbative, in the sense that the configuration of the non¬ 
linear held at the end of the process is very different from 
what it was at the beginning. Here we have in mind the “stan¬ 
dard” perturbation theory, rather than the one based on the 
Inverse Scattering Transform ll2^ (indeed, below we will dis¬ 
cuss and use the latter). In our numerical experiments, we 
initialize a “breather” state ll24l - lz^ (see Eig. [T]), which is a 
multi-soliton solution of the NLSE where the centers and the 
velocities of the constituent solitons coincide at the initial time 
(unlike in the more widely-known case of the breather in the 
sine-Gordon equation, where the breather is a distinct kind of 
solution from the two-soliton solutions IJt]). We then send it 
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FIG. 1. The 1:3 odd-number ratio breather [see Eq. 0], at three points 
of time. Here T = 327r (in natural units; see text) is the period of density 
oscillations. This is a particular two-soliton solution of the NLSE, where the 
constituent solitons have coinciding centers and zero velocities; the norms of the 
constituent solitons are 1/4 and 3/4. 


to collide with a barrier—a repulsive potential bump, in our 
case either Gaussian or rectangular. As Fig. |2] shows, there 
exists a parameter regime in which the scenario outlined in 
the introductory paragraph takes place. In all cases, the sizes 
of the emerging solitons are precisely the sizes of the soliton 
constituents of the original breather, so that the transmition 
of each individual soliton is all or nothing. This results in a 
staircase-like transmission plot of Fig.|3] The cases where at 
least one soliton is fully reflected are particularly clear indi¬ 
cators that the process is nonperturbative (again, in the “stan¬ 
dard” sense of that word). The decomposition of breathers 
has been studied previously, but, at least when it comes to 
decomposition into the constituent 1-solitons, only for weak 
integrability-breaking perturbations 12512^ . 

Details of our system are as follows. In the context of ultra- 
cold Bose gases, the ID nonlinear Schrodinger equation (also 
called the Gross-Pitaevskii equation in that setting) is 

ihdt^ = giDiV,|T'pT' + y(z)T-, (1) 

2m 

where \^{z,f)\^dz = 1. It describes the dynam¬ 
ics of the order parameter f) of the Bose-Einstein 
(quasi)condensate of A), atoms of mass m confined to a cigar¬ 
shaped waveguide, in which the transverse confinement (a 
harmonic potential, imw^r^) is so tight that the energy spac¬ 
ing hwr between the ground and the first excited transverse 
state is much greater than the typical energy of the confined 
particles, making the system effectively ID. The atoms inter¬ 
act with an effective pairwise interaction giu 5{zj — Zk), where 
Zj and Zk are the positions of the atoms. For us, gio < 0. 
Ns, is assumed large enou^ that we expect that quantum cor¬ 
rections to be negligible 1^ (we thus ignore the fact that 
it really should enter as — 1). The potential V(x) will 

be our barrier: V{z) = Vb exp(—2(z/«))^), where w is the 


a £kin/V6 =0.296875 



b £kin/l« =0.30625 




FIG. 2. Three kinds of outcome following a collision of a 1:3 ONRB with 
a Gaussian barrier. All numbers are in natural units (see text). The barrier is 
at ir = 0, has peak height Vo, and 1/e^ half-width w = 0.9. The breather 
has impact speed of 0.025, and chemical potential which is much higher than its 
(center-of-mass) kinetic energy per particle. The main plots show the details of the 
collision process, where the various density profiles correspond to different times 
£; the insets are the final outcomes, a, total reflection; b, the constituent soliton 
of norm 1/4 is transmitted, while the soliton of norm 3/4 is reflected; c, both 
constituent solitons are transmitted. The underlying experimental values are as in 
the current lithium experiments, except for the impact velocity and barrier width, 
which are ten times lower and larger, respectively. This boosts the ratio of the 
chemical potential to kinetic energy per particle, without which the conservation 
of the shapes of the individual solitons would not be as perfect (see Extended Data 
Figure 1 in the Supplementary Information). In the flnal outcomes, the constituent 
solitons are well-separated even when they are on the same side of the barrier 
because they generally spend unequal times on top of the barrier, and also because 
they generally leave the barrier with unequal speeds. 










































half-width of the beam at the 1/e^ intensity point. We have 
numerically verified that a rectangular barrier of comparable 
dimensions gives qualitatively the same results. If V(z) = 0, 
the equation is integrable, both on the whole real line and on a 
ring (i.e. with periodic boundary conditions) l[ll]. Best known 
are the exact solutions on the real line, which may be obtained 
through the inverse scattering transform [[llllSl^ which we 
will use in our work. It is convenient to work in the “natural 
units” (see the Supplemetary materials) in which the equation 
takes the form 

idttp = -^d^ip-lipl^ip + V(z)ip, ( 2 ) 


t)\'^ dz = Af, where the new normalization J\f 
may be chosen at will. Our parameters are chosen so that 
the system is in a regime that is, on the one hand, as close 
as possible to that obtained in actual experiments with lithium 
11221] , but in which, on the other hand, we do see the “staircase” 
transmission plot. Below we will discuss how to observe this 
experimentally. 

Let us consider the integrable case, V{z) = 0. The ba¬ 
sic 1-soliton solution, of norm 1, is ^ exp(i</8) sech( 0 / 2 ). 
Note that if t) (of norm Af) is a solution of Eq. (O, 
then so are the Galilean-transformed solution e'yq>\iv{z — 
zq) — ip{z — vt — zq, t) and the renormalized solution 

of norm ^Af', these two transformations may be 
used to produce 1-solitons of arbitrary norm, initial position, 
and velocity. iV-soliton solutions (see the Supplementary Dis¬ 
cussion for details on what follows) are parametrized by AN 
real parameters, four per constituent soliton Aj, vj, Zj^, 
and j = 1, .. ., N. The norm of ijj is at 

t —>■ ±oo, if the velocities vj are all distinct, ip becomes a sum 
of N 1-solitons, of norms Aj, traveling at velocities Vj', Zj^, 
and 0 almost the position and phase, at f = 0, of the 
jth soliton. We do not consider degenerate cases when two 
or more constituent solitons have both the same norm and the 
same velocity. 

In the numerical experiments discussed so far, the initial 
state of ■!/; is a (Galilei-boosted) two-soliton breather belong¬ 
ing to a particular sequence multisoliton solutions which we 
will call the odd-norm-ratio breathers (ONRBs) j^ . The nth 
element of this sequence consists of n overlapping solitons, 
of norm ratios 1:3: • • • : (2n — 1), which periodically 
beat (“breathe”) in space due to interference and the fact that 
the angular velocities of the phases of individual solitons de¬ 
pend on their norm. All the ONRBs have the property that at 
least once per breathing period, their waveform is a sech, i.e. 
the same as a one-soliton solution for a NLSE with a different 
coupling constant. As we will discuss below, this provides a 
way to excite them experimentally. The two-soliton ONRB 
has the form 


■*/'ONRB[2](a:, t) 


e® 128 ^cosh ^ -I- 3 e® ^ cosh 
6 cos + 8 cosh f -I- 2 cosh | 


(3) 


it may be obtained from the general two-soliton solution dis¬ 
cussed above if one sets Ai = 3/4, A 2 = 1/4, i/i o = 0, 
(f>2,o = TJ'j and vi = V2 = 2:1,0 = 2:2.0 = 0. The breather 



FIG. 3. Figure 3 | Transmission plot for the scattering of a 1:3 ONRB 
off a Gaussian potential. The dash-dotted horizontal lines are at transmissions 
values of 1/4 and 1, corresponding, respectively, to only the smaller soliton being 
transmitted, and to both solitons being transmitted. The three curves correvSpond to 
three different values of barrier width w; in all cases, as a function of -Ekin/Vo, the 
transmission proceeds in the same three steps a, b, and c presented in Fig.[^ with a 
substantially wide plateau during which transmission is 1/4. All other parameters 
are as in Fig.[^ See Extended Data Figure 4 in the Supplementary Information for 
how this plot would look if the constituent solitons were completely decoupled. 


is sometimes referred to as a bound state of the two solitons 
ll24l] . but note that the dissociation energy is zero 1^ : 
if, instead of vi = V 2 = 0, we have jui — t; 2 | = e > 0, 
no matter how small, then for t ^ 1/e, the solitons will be 
well-separated and ballistically receding from each other with 
a speed that tends to e as f —oo: constituent solitons of true 
multisoliton solutions experience no mutual force. (This is 
in contrast to the well-known force between simple soli¬ 
tons, whose sign depends on relative phase; but in that case, 
the initial state is not an exact multisoliton solution but rather a 
linear superposition of simple solitons, relevant when simple 
solitons are injected into the system from the outside.) The 
zero dissociation energy also means that the internal kinetic 
and internal potential energy are fine-tuned, due to integrabil- 
ity, to exactly cancel. We have also made sure the breather 
always impacts the barrier at the same phase in its breathing 
cycle. 

Several additional points are in order: 

(i) We have checked that once the solitons move away from 
the barrier, the field ijj is numerically of the solitonic form; 
e.g., when the constituent solitons are spatially separated, their 
form is that of 1-soliton solutions, Galilei-boosted, of norms 
1/4 and 3/4. 

(ii) From Fig. |2] we see that the nonlinear interactions are 
not negligible at any point during the collision, since the field 
density is not close to zero anywhere over the barrier. 

(iii) The stable plateaus in the transmission plots disappear 
if the underlying equation is not integrable; see Extended Data 
Figure 1 in the Supplementary Information. 

(iv) The phase in the breathing cycle at which the collision 
occurs does affect where the plateau at 1/4 begins and ends, 
but not its existence; see Extended Data Figure 2. 
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(v) Similar results may be obtained with multisolitons of 
any number of constitutent solitons; for example, superposi¬ 
tions of more than two solitons also transmit constituent soli¬ 
tons “one by one,” resulting in multiple plateaus; see Extended 
Data Figure 3. 

Here is why the effect happens. Given the NLSE of Eq. (|2|l, 
the inverse scattering transform (1ST) introduces, for each 
time t, an auxiliary linear problem, in which the field '0(z, t) 
of the NLSE plays the role of a potential: 

f -id^ 'tp*{z, t) \ f u(z, t) \ _ u{z, t) \ 

-il^iz, t) id^ J v{z, t) ) u(z, t) J 

(4) 

The matrix on the left-hand side is called a Lax operator; 
the scattering data for this problem consist of: 1. the so- 
called reflection coefficient; 2. the normalization coefficients 
(which are really the prefactors of the long-distance asymp¬ 
totics of the bound states); and 3. the discrete eigenvalues A, 
which are the only component of the scattering data that will 
enter our analysis below. Scattering problems of this kind can 
formally be introduced for any nonlinear PDE, but only for 
an integrable PDE can the reflection and normalization coef¬ 
ficients of the corresponding Lax operator be propagated in 
time easily (using a second Lax operator, which we will not 
discuss further), while the discrete eigenvalues A remain time 
independent O^ . In fact, it turns out that each distinct A de¬ 
scribes a constituent soliton of the held ip{z, t): up to numer¬ 
ical prefactors of order 1 in the natural units, the real part is 
the soliton velocity v, and the imaginary part, the norm A (see 
above). The reflection coefficient, in turn, corresponds to the 
radiation component of the held ip, and is strictly zero if ip 
consists purely of solitons. To solve the initial value problem 
for the NLSE, one proceeds in a way that has often been called 
a nonlinear generalization of solving a PDE via Fourier trans¬ 
forms, where the scattering data play the role of the Fourier 
coefficients: one first finds the scattering data corresponding 
to the “scattering potential” ip{z, t = 0); this is a nonlinear 
analog of taking the direct Fourier transform. One then (eas¬ 
ily) propagates the scattering data forward to time t. Finally, 
one solves the inverse scattering problem for the propagated 
scattering data to obtain the “scattering potential” ip{z, t) at 
any later time t. This final part reduces to a linear integral 
equation, which is analytically tractable in many physically 
important cases, and which is the nonlinear analog of taking 
the inverse Fourier transform. 

In the presence of an integrability-breaking potential 
V{z) = evext.i^), where VextXz) ^ 1 and e positive but not 
(yet) assumed small, the As develop time dependence, which 
is described by certain well-known relations of so-called 1ST 
perturbation theory ll^ . For our purposes, however, rather 
than using these relations in the form in which they are usu¬ 
ally presented, it is more convenient to base our analysis on 
the following exact relationship, for which we provide, in the 
Supplementary Discussion, an alternative derivation based on 
the Hellmann-Feynman theorem: 

dX _ .e dzve^t. {u"^ Ip + Ip*) 

^ “ *2 r dzuv ’ ^ ^ 

J —OO 


where ip{z, t) is the exact solution of the perturbed NLSE. 
The main difference as compared to the usual form of this 
relation S is that instead of bound-state asymptotics one 
has a straightforward normalization integral. The total change 
in A is A A = (dX/dt) dt. 

In this language, the effect we are studying—the preser¬ 
vation of soliton identity through the collision—translates to 
the imaginary parts of As remaining unaffected by the col¬ 
lision, while their real parts (soliton velocities) change sub¬ 
stantially. To see the mechanism for this, let us estimate how 
AA scales with e In the regime of interest, the norm of the 
breather (and so of its constituent solitons) is ~ 1, while the 
impact kinetic energy of the breather is of the same order as 
the potential, meaning that the impact velocity of the breather 
~ <C 1. Note that the barrier appreciably affects the so¬ 
lution only while the breather is on top of it. To proceed, 
we make the following assumption, which seems very phys¬ 
ically reasonable even though we cannot, at present, justify 
it with any degree of rigor: as long as the sizes of the con¬ 
stituent solitons (the imaginary parts of the discrete eigenval¬ 
ues of the Lax operator) are not much different from their 
initial values, the centers of mass of the constituent solitons 
will, at least qualitatively, move as if they were non-mutually- 
interacting particles undergoing a classical collision with an 
effective barrier (a convolution ll2^ of the actual barrier and 
\ip\'^), even while the constituent solitons are not well sepa¬ 
rated. The typical collision time of a single soliton can be 
estimated as follows: since the barrier is much narrower than 
the breather, the width of the effective potential is the same as 
that of the breather, i.e. ^ 1. Thus the time that constituent 
solitons spend atop the barrier scales as width/speed ^1/ y/e. 
But then AA ^ i/e. Now recall that A ^ velocity -f inorm; 
thus, AA is comparable in magnitude to the velocity, but is 
much smaller than the norm, which thus remains unaffected, 
for all intents and purposes. And in the case when one of the 
solitons has just the critical velocity to spend a long time atop 
the barrier, the other solitons, having different norms, will not, 
and will leave the barrier in the time we just estimated. Thus, 
as the breather begins to climb the barrier, the changes in the 
norms of the constituent solitons (imaginary parts of their As) 
begin to accumulate according to Eq. Q; but the solitons will 
leave the barrier before the change in their norms can become 
significant—with a possible exception of one of the solitons 
having just the right size-velocity combination to linger atop 
the barrier. However, since the other solitons leave the barrier 
quickly, the remaining soliton is energetically protected from 
fragmenting. 

In the usual sense of the word, the process is heavily non- 
perturbative, since it results in a dramatic change in the shape 
of the held ip. Almost all of the parameters that determine 
the properties of the constituent solitons substantially change 
as well, save for the norms. If we metaphorically picture the 
perturbation as a kind of heating up of the previously con¬ 
served (“frozen”) quantities, then we may say that the norms 
are “superheated”: though one might think they should start 
flowing, they do not. And the reason is that the 1ST, counter¬ 
intuitively, treats the soliton velocity and norm as a unit, as the 
real and imaginary parts of a coordinate A in the “IST space”: 
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the magnitude of a change in A may be much smaller than 
the imaginary part of A, while being substantial compared to 
the real part. The “superheating” of the norms is thus a direct 
manifestation of the underlying 1ST structure of the problem. 

Let us discuss the possible experimental realization of the 
effect we have just described (we thank R. G. Hulet’s group 
for many useful discussions and for providing us with their 
current experimental parameters). The protocol for the cre¬ 
ation of the 1:3 breather is based on the observation that, at 
f = 0, its wavefunction coincides with the wavefunction of 
the stationary soliton but for the NLSE whose coupling con¬ 
stant is one quarter the coupling constant for the actual NLSE 
mi: ^/’oNRB[ 2 ](a^, 0) = isech (|). Since a single soliton is 
usually first created very close to the stability threshold 
and the breather is even more unstable than a single soliton, 
the magnitude of the coupling constant should first be adia- 
batically reduced by a factor of ten, and only then quenched 
up by a factor of four (rough variational estimates suggest that 
the resulting breather should be well within stable regime). In 
the ^Li experiments 123] . the parameters of Eq. ( IS2l i have the 
following values: = 27r x 254 Hz; ojz = 27r x 31 Hz; 

giD = 2HuJra, where a is the (3D) scattering length of the 
atoms, experimentally —1.0 x gq (uq being the Bohr radius). 
It can be tuned down in magnitude by at least a factor of 10 by 
using the magnetic Eeshbach resonance = 3 x lO'^, 

which is indeed at the collapse threshold. The barrier is pro¬ 
duced by a 900 GHz blue-detuned Gaussian laser beam result¬ 
ing in the half-width w of the beam at the 1 /e^ intensity point 
of re = 4.2 fim, while t^//i is of the order of 200 — 3000 Hz 
(where h is the Planck’s constant). In experiments, V(z) will, 
in addition to the barrier, also include a longitudinal harmonic 
confinement we will omit it from the analysis, but 

we have checked that if its value what it in the experiments 
(see below), it does not qualitatively affect our results. The 
typical impact speed in the experiments is 514 pm/s. 

In our numerical experiments, we took a = — 0.41ao; in 
order to produce very sharp transitions in transmission plots, 
we lowered the impact velocity tenfold from the typical ex¬ 
perimental value of 514 /rm/s. Eor the smaller soliton, this re¬ 
sulted in (f?rest/.^a)/E'K.i ~ 8. However, stable plateaus can 
be obtained even with larger velocities, though the transitions 
will no longer be as sharp. Also, to prevent the computational 
grid from becoming impractically big, we increased the width 
of the barrier tenfold from the experimental value (results of 
Fig. |3] suggest this will not matter qualitatively). 

When the superheating of soliton norms does become 
experimentally accessible, the existence of a transmission 
plateau could be used to stabilize bright solitonic atom inter¬ 
ferometers against the fluctuations of the profile of the beam 
splitter potential, and against the parasitic mean-field effects 
1^1^. In the latter case, the problem is that due to interatomic 
interactions, the amount of phase that a matter wave accumu¬ 
lates as it travels depends also on its size. But the current 
interferometric schemes rely on splitting a single soliton off 
a barrier to send a portion of it down one interferometer arm, 
and a portion down the other. So a lack of control over what 
portion of atoms goes to each arm results in the accumula¬ 
tion of different amounts of phases (which one does not know 


how to account for) in each arm during free flight jst]. For 
flight times that will be necessary to produce competitive ac¬ 
celerometers, the lack of control present in current schemes 
will be unacceptably large. An interferometric scheme that 
takes advantage of the superheating of soliton norms, how¬ 
ever, may be able to overcome this problem. 

Let us conclude with an outlook as far as superheated inte- 
grability (SI) as such. Its essential ingredients are (i) a sepa¬ 
ration of scales in the magnitudes of the real and the imagi¬ 
nary parts of the eigenvalues of the Lax operator; and (ii) an 
integrability-breaking perturbation whose effect is, for each 
eigenvalue A, greater than or comparable to the smaller of the 
pair (|ReA|, |ImA|), but much smaller than the greater of that 
pair. Let us call this the ’’superheating regime” of the per¬ 
turbation; the nontrivial aspect of predicting or explaining an 
occurrence of SI is estimating the magnitude of the effect of 
the perturbation. But once these essential ingredients are un¬ 
derstood, it becomes clear that there must exist other types 
of SI, in other settings, even including other underlying inte- 
grable systems. In fact, we can already give another example 
of SI in the NLSE, in which integrability is broken not by a 
collision with a barrier, but by phase imprinting; see Extended 
Data Figure 5 in the Supplementary Information. 

We thank B. Sundaram, B. A. Malomed, R. Hulet, and R 
Dyke for useful discussions; we additionally thank R. Hulet 
and P. Dyke for providing us with detailed parameters of their 
experiments. This work was supported by grants from the Of¬ 
fice of Naval Research (N00014-12-1-0400) and the National 
Science Foundation (PHY-1402249). 


Supplementary Discussion 


1. Energetically allowed eractional loss oe par¬ 
ticles EOR A SINGLE SOLITON 

Consider the stationary, single-soliton solution (normalized 
to 1) of Eq. (1) with U = 0 in the main text, 'I'i(z) = 
y^c/2sech(cz), where c = |pid| AfaTO/(2/i^). 

The energy-per-particle of the soliton is given by 

E/N., = dz + 

= -dNl withd=^, 

where pio < 0. A fracturing of a single soliton into two 
smaller ones, Ni + N 2 , will raise the energy-per- 

particle by 2dNiN2. If there are multiple fragments, their en¬ 
ergies are even less negative, so the increase in the energy-per- 
particle will be even greater; same if some of the particles are 
radiated away, as their energy will actually be positive. So the 
assumption of a fracture into exactly two smaller solitons will 
give the maximal possible particle loss. If the fracturing is to 
occur due to the soliton (with center-of-mass speed v) collid¬ 
ing with a barrier, this energy can come only from the kinetic 
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energy per particle of the center of mass, i?K,i = ^mv^. Let 
us write Ni = N^ — AN and N 2 = AN, where, without loss 
of generality, we impose the condition that 0 ^ AN < N^/2. 
We therefore have 0 ^ 2d(A^a — AN)AN ^ dividing 
through by 2dN^, we get 


where the Uj’s are the solutions of the system of N equations 


E 


+7fc 

^3 + K 


Uk = l, 


j = 1, 


N, 




AN 


l;k,i 

nJ " 2dN^ • 


(SI) 


where 


^3 — /2 + iuj 


Now we invoke the assumption that E^k,! 'C |E/iVa|. Thus 
the right-hand side of the above relation is much less than one. 

The equation x—x^ = e has the solutions xi^2 = ~ 

i.e. one close to zero, the other close to one. Because we 
imposed that AN < Nn/2, it follows that AN/N^ must be 
much less than one; so we may neglect the quadratic term in 
Eq. (ED, obtaining 


AN 


-E'k,! 
2dN^ ■ 


See also Ref. (the paragraph before Sec. 3). 


and 


= exp [Xj{z — Zjfi) + iA|f/2 + i(l>j,o] ■ 

Here the parameters zj^q, and (f)j Q are almost but not quite 
the position and phase, at f = 0, of the jth soliton when it 
is spatially separated from the others (see below), but the first 
two are precisely its norm and velocity, as we now explain. 
The norm of ijj is as f —>■ ±c», if the velocities 

Vj are all distinct, tp becomes a sum of N 1-soliton solution, 
of norms Aj, traveling at velocities vj. More precisely, in 
the limit as the jth soliton becomes more and more spatially 
separated from the others, its form converges to 


2. The natural units for the ID nonlinear 

SCHRODINGER EQUATION 


-A sech 



Z3) + q, 


exp[i(0j -f T-j)], 


Consider the ID nonlinear Schrodinger equation 

ihdt'H = + PmiValT-l V + V{z)'i>, (S2) 

Zm 

where t)\'^ dz = 1. By choosing to work in “nat¬ 

ural units,” it may be written as 

idtip =± \ip\^p) + V{z)^p, 

where the sign of the |'0p'0 term is the sign of pio. In the 
present work, this is negative (i.e. attractive, or “focusing”), 
so that it can support bright solitons. ^ = 

J\f is the new normalization, which may be chosen at will. 
In the resulting system of units, the unit of length is ul = 
Afb?/{mgiuNn), of time, ut = /(mgij^N^), and the 

field is transformed as T'(z, t) = ’p^zjuL-, tjuT)!\/JT ul- 


3. The A-soliton solution of the nonlinear 
Schrodinger equation 

The multisoliton solutions of the equation 

idtipiz, t) = -]^dlip{z, t) - \3l3{z, t) 

may be constructed as follows IHt]: given the 4iV real param¬ 
eters Aj > 0, Vj, Zjfi, and <j)jfi, j = 1, ■ • ■, -A, we have that 

N 

Hz, t) = J2uj{z, t), 

3 = 1 


where 


Zj = Zj,o + Vjt, 

H = ^3iz - Zj) + i (^ 2/4 -f v]) t + , 

and the real numbers qj and rkj capture the interaction with 
the other solitons. They are given as 


qj + mij 


N 

^sign(zfc - Zj)\n 


k=l 




Aj Ak 2,%{vj 
Aj — Ak -\- 2i(vj 


Vk) 
Vk) ’ 


where signz is -1, 0, or -i-l if z < 0, z = 0, or z > 0, 
respectively. 

Now we see that the parameters Zj^, and pj^o differ, 
through qj and T'j, respectively, from the actual position and 
phase at f = 0 of the jth soliton when it is spatially separated 
from the others. But the actual position and phase at f = 0 is 
easily characterized; since the displacements qj depend only 
on the ordering of the spatial positions of the solitons (and not 
on the magnitudes of the relative distances), it follows that if 
one wants isolated solitons to sit at zqj at f = 0, one may 
proceed as follows: first set each Zjfi to zoj, and compute the 
gj’s. Then set Zj^o = zqj + 2qj j Aj, and proceed to compute 
the A-soliton solution; a similar correction may be applied for 
the initial phases. 

The case of degenerate \j, when two or more constituent 
solitons have both the same norm and the same velocity, may 
be treated by taking the appropriate limit of the system of N 
equations above. In this case one obtains solutions that are 
qualitatively different from those discussed thus far. For ex¬ 
ample, in the two-soliton case, as f ± 00 , one finds 
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that the distance between the solitons increases proportionally 
to hi{A^t) (in natural units). Since here the solitons separate 
on their own on the time scale of \jA?, the collision exper¬ 
iment should last shorter than that. On the other hand, the 
breather needs to start sufficiently far from the barrier so that 
it begins in an approximately integrable regime, and it needs 
to be sufficiently slow so that the kinetic energy per particle is 
much less than the chemical potential. It turns out that these 
constraints are impossible to satisfy simultaneously, and thus 
the degenerate case is not of interest for us. 


4. Derivation of the exact expression for dX/dt, 
EQ. (5) IN THE MAIN TEXT 


where 



To convert expressions in this text to the ones appearing in 
Ref. 1^ one should use the following replacement table; 


z —>■ x 

'>p{z, t) —>■ u{x, t) 
u{z, t) —>• t) 

v{z, t) —>• 'lp^‘^\x, t) 

Vext{z, t)tl;{x, t) -)> i{P[u]){x, t). 


We are dealing with the ID nonlinear Schrodinger equation 
in Eq. (1) in the main text. 


d h? 

+ siDW.|a'(2, t) + v{2, t), (S3) 

where pio < 0, in the presence of the external integrability- 
breaking potential V{z, t) (note that here we will allow this 
external potential to explicitly dep end on time). In order to 
facilitate comparison with Ref.|^ which treats the same kind 
of problem, we will work in the units in which h = 1, m = 
1/2, and |piD|.A^a = 2. Thus the NLSE becomes 


d . . . 

■s'*'"- = 


o2 

OP 0(^,0 


+ eVext{z,t)'fp{z, t) , (S4) 

where the external potential in (IS31 l is factorized as 

V(z, t) = VoVext{z,t) 
max[vext{,z, 0)] = 1, 

Z 

and the small parameter e is 

^ 2h'^Vo 

The hrst Lax operator reads 


£ = 


L M 
-Mt -L 


with 


f 

M = ^*{z,t) . 


(55) 

(56) 

(57) 


Eor each instance of time t, one can set up an eigenstate- 
eigenvalue problem, which is the central object of interest of 
this derivation; 


4.1. Relevant functional analysis 


Erom the fact that L is Hermitian, U = L, it follows that 
the Lax operator (IS51 l possesses the following property; 

. 

This property induces a particular Hermitian form (•, •), a 
pseudo-inner product; 


\W2>-) =<Wi\w2)^= {Ui\u2) - (Ui|t;2) 

= Jdz {ul{z)u 2 iz) - vl{z)v 2 iz)} . (S8) 

This Hermitian form lacks the property of being positive def¬ 
inite (i.e. lacks the property that, for all \w), (ti;|ti;) > 0 and 
(tultu) = 0 if and only if jw) = |0)). The rest of the inner 
product axioms, on the other hand, remain intact; 

-<W2\wiy = ^W2\wiy* 

and 

-<wi\aw2 + bw3)^= a^wi\w2>- +b^wi\w3>- . 

The Lax operator £ from Eq. (IS5I) is symmetric with respect 
to this form; 

£|w2 )-) = (£|wi:^, . (S9) 

The property above justifies a standard notation 
C\w 2 >-^ =< wi\C\w 2 >■ that we are going to 

employ below. 

The pseudo-Hermiticity property (IS9I ) implies the follow¬ 
ing properties of the eigenstates of £; Let C\wi)^= Ai |wi;^, 
C\w 2 >-= A 2 |w 2 :^, and C\w>-= A|i(;;^. Then 

1 . eigenvectors whose eigenvalues are not complex conju¬ 
gates of each other are mutually orthogonal. 


t\w)^= Alw;^ 


•^1 7 ^ ^^2 ^ -<Wl\w2i^= 0 
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2 . non-zero norm eigenstates of L correspond to real 
eigenvalues (a corollary of the above); 

A* A => ^w|w:^=0. 


Note that the eigenspectrum of L is not necessarily complete. 
In cases when the Lax operator (IS51 l represents a linear stabil¬ 
ity analysis equation of a nonlinear PDE, the missing states are 
associated with the continuous symmetries of the PDE that is 
broken by the solution in question ifs^l . (Eor a “flat” conden¬ 
sate, = const, we found one missing state; there could 
be more. There seem to be none for a single soliton.) 

The operator ( IS5b also possesses properties specific to a 
particular form of the matrix elements, Eqs. (IS61 l and (IS71 i: 

L* = -L . 

This property implies that; 

1. real eigenvalues A are doubly degenerate. The corre¬ 
sponding eigenstates. 


|w:^= 


|w:^= 


u{z) \ 
v(z) ) 

u{z) \ 
viz) ) 


are related by 


u{z) = —v*{z) 
v{z) = +u*{z) . 


Within the context of the Inverse Scattering Transform, the 
wavefuction t) in the parent NLSE, Eq. (IS3b , is assumed 
to be localized in space, while the eigenstates of the Lax op¬ 
erator £ of Eq. (IS5b are required to be finite at a; = ±oo. In 
this case, the real eigenvalues of £ form a continuum spec¬ 
trum, while the complex eigenvalues are discrete. Einally, in 
the parent NLSE, the complex eigenvalues correspond to the 
solitonic part of the scattering data, while the real eigenvalues 
correspond to the thermal noise. 

Erom now on, we will assume that the eigenstates of £ with 
substantially complex eigenvalues (i.e. the “discrete spec¬ 
trum”, or “bound states”) are normalized as 


-<w-\w+y=l. (Sll) 


Eor the case of a single soliton, 

ip{z) = —i sech(a::), 


the corresponding eigenvalues and eigenstates are 


and 



_t_2 

|w+^= — exp(a;/2) 


-1-1- tanh(a:) 
sech(a:) 



— 2 

|u;_^= — exp(a;/2) 


—sech(a;) 
-1-1- tanh(a;) 


Here, 

C\wy= A|w:^ 

C\wy= A|w:^ . 

2. Complex eigenvalues A come in complex conjugate 
pairs, A+, A_ such that A_ = (A+)*. The correspond¬ 
ing eigenstates. 



are related by 

u-{z) = -{v+y{z) 

v-iz) =+iu+nz) . (SIO) 


4.2. The Hellmann-Feynman theorem 

Let £ depend on a parameter £ = £(^). Its discrete 
eigenvalues A± and the corresponding eigenstates, then 
also depend on A± = A±(^), and |w±(^);^. 

Let us express the eigenvalue A+ as 

^+(0 =<w-{Cj\C{£,)\w+{y)>- . 

Then, the derivative of A+ with respect to ^ becomes 

= (I'”’-'-- 

+ [fy] '"’+=■)+ 

As in the proof of the usual Hellmann-Eeynman theorem, the 
sum of the first and the last term will turn out to be propor¬ 
tional to the derivative of the norm; and since the norm is one, 
the derivative is zero. Indeed, the first term gives 


A+|w+^ 



Here, 


£|rt;+:^= A+|rt;+:^ 
£|rt;_:^= X-\w->- 
A_ = (A+)*. 
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Similarly, the last term gives 


w->-, C — \w+>- ) = C\w->-, —\w+>- 


dC 


= ( A_|u>_^, 


dC 


= (A-)* —\w+'^] . 

But (A_)* = A 4 .; thus, the sum of the two terms gives 


-|»+> 


Notice that the contribution to dX/dt from F[tlj] disappears. 
Indeed this contribution describes the time derivative of the 
Lax eigenvalue in the time evolution according to the unper- 
turbed^FS', this derivative indeed vanishes as a consequence 
of integrability of the NLSE. 

A translation to the Kivshar-Malomed notation system of 
Ref. I 23 gives 


^A _ 1 1 


dz t, Xn)eP[tfj]{x, t) 

t)e*P*[i)\{x, f)| . 


Thus, we get the following generalization of the Hellmann- 
Feynman theorem: 



Au>_| 


(|t)!«.+!- 


(S12) 


4.3. The exact expression for dX/dt from the Hellmann- 
Feynman theorem 


Let us set 

m = 


-i-§^ t) 

-tpiz, t) +z^ 

C{t)\w{t)>-= X{t)\w{t)>- 


dt 


m = 

-{F[f\+eP[f]){z,t) 0 

f) 


v{z, t) 

/ -\-oo 

u{z, t)v{z, f) = 1 , 

-00 


where F[tf}]{z, t) = —i —-^ — 2\ip{z,t)\‘^ tf}{z, t), and 

P['ip]{z, f) = —ivext{z, t)'>P{z, t). According to the 
Hellmann-Feynman theorem, Fq. (IS12b . the time derivative 
of the Lax eigenvalue is 




x = 


^ + CXD 


- dz {uliz, t)eP[ip]{z, t) 

J —00 

-vl{z, t)e*P*[f]{z, <)} 

/ +00 

dz {eu^(z, t)ip(z, f) 

-00 

+e*t;^(z, t)ip*{z, t)}vext{z, t). 
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Extended Data 



Extended Data Figure 1. Integrable vs. nonintegrable case. The initial state for the nonintegrable case was prepared by time- 
propagating the breather at rest while the nonlinearity was slowly ramped from |^/)p to with p = 3/2. The result was Galilei- 
boosted and scattered off of a Gaussian barrier whose width (for numerical reasons) was twice the reference value. All other parame¬ 
ters were at their reference values, including the number of particles. Also plotted is the integrable case, all of whose parameters are 
at their reference values. In particular, the number of particles is four times smaller than that for Fig. 3 in the main text, resulting in 
a degraded, but still noticeable plateau at around 25% transmission. In the nonintegrable case, no plateau can be discerned. 







transmission 



Extended Data Figure 3. Transmission plot for a three-soliton breather solution. The constituent solitons have norms 1/6, 1/3, and 1/2 
(1:2:3 norm ratio, which does not belong to the sequence of odd number ratios). The dash-dotted horizontal lines are at transmissions 
values of 1/6, 1/2, and 1, corresponding, respectively, to only the smallest, norm-1/6 soliton being transmitted, to the norm-1/6 and 
norm-1/3 solitons being transmitted, and to ail three constituent solitons being transmitted. 
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0.2 0.4 0.6 0.8 1.0 1.2 1.4 



Extended Data Figure 4. Uncoupled solitons vs. breather. How the “staircase” plot of Fig. 3 in the main text would look if the 
constituent solitons of the breather were completely uncoupled, as compared to now it is in reality, a, The transmission plots for the 
scattering of single solitons, of norms 1/4 and 3/4, off a harrier. All parameters are as in Fig. 3 in the main text, with harrier width 
w = Wo. h, dashed line: the weighted sum of the single-soliton transmission curves from panel a, with the norms used as weights. 
Solid line: the transmission curve for the breather for the same set of parameters. This is the same curve as the w = wo curve in Fig. 3 
in the main text. 
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Extended Data Figure 5. Superheated integrability in the context of phase imprinting. At t = 0, the breather wavefunction is multi¬ 
plied by the space-dependent pure phase Here e = 0.01 andv5(a:) = [c- cos(2Tvmx/L) + Sm sm(27rma;/L)], 

with 1/ = 16 and M = 5; Cm and Sm were drawn from the uniform distribution on [—1, 1], and in the realization shown here had 
the values (ci, ..., cs) = (0.307, 0.622, 0.648, -0.738, 0.304) and (si, ..., ss) = (0.353, -0.0422, -0.794, -0.746, 0.721). The 
function is plotted in the inset. The main plot shows the time evolution of the density |^/)|^, during which the constituent 
solitons separate. Once they are well-separated, one can verify that their norms are 0.25 and 0.75. In the context of Eq. (5) 
in the main text: if one uses the approximation « 1 + ieip{x), then the process depicted in this Figure corresponds to 

Vext.{x, t) '4){x, t) = i5{t) tp{x) tp{x, r). Just as in the case of collision with a barrier, the separation of scales between the 

real and Imaginary parts of the Lax-operator eigenvalues A, which correspond to the constituent solitons, again imply that the sollton 


velocities Re A) change but norms (~ Im A) do not. 
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